A Quantitative Data-Driven Analysis Framework for Resting-State Functional Magnetic Resonance Imaging: A Study of the Impact of Adult Age

The objective of this study is to introduce a new quantitative data-driven analysis (QDA) framework for the analysis of resting-state fMRI (R-fMRI) and use it to investigate the effect of adult age on resting-state functional connectivity (RFC). Whole-brain R-fMRI measurements were conducted on a 3T clinical MRI scanner in 227 healthy adult volunteers (N = 227, aged 18–76 years old, male/female = 99/128). With the proposed QDA framework we derived two types of voxel-wise RFC metrics: the connectivity strength index and connectivity density index utilizing the convolutions of the cross-correlation histogram with different kernels. Furthermore, we assessed the negative and positive portions of these metrics separately. With the QDA framework we found age-related declines of RFC metrics in the superior and middle frontal gyri, posterior cingulate cortex (PCC), right insula and inferior parietal lobule of the default mode network (DMN), which resembles previously reported results using other types of RFC data processing methods. Importantly, our new findings complement previously undocumented results in the following aspects: (1) the PCC and right insula are anti-correlated and tend to manifest simultaneously declines of both the negative and positive connectivity strength with subjects’ age; (2) separate assessment of the negative and positive RFC metrics provides enhanced sensitivity to the aging effect; and (3) the sensorimotor network depicts enhanced negative connectivity strength with the adult age. The proposed QDA framework can produce threshold-free and voxel-wise RFC metrics from R-fMRI data. The detected adult age effect is largely consistent with previously reported studies using different R-fMRI analysis approaches. Moreover, the separate assessment of the negative and positive contributions to the RFC metrics can enhance the RFC sensitivity and clarify some of the mixed results in the literature regarding to the DMN and sensorimotor network involvement in adult aging.


INTRODUCTION
Among the different analysis approaches for resting-state fMRI (R-fMRI) data, the anatomic region-of-interest (ROI)-based, and data-driven independent component analysis (ICA) methods are probably the most commonly used (Tyszka et al., 2014). Restingstate functional connectivity (RFC) results from the ROI-based and ICA derived methods are generally similar but conceptually different. The quantitative relationship between ROI-based and ICA derived measures of RFC has been investigated with computer simulation and experiment approaches (Joel et al., 2011;Rosazza et al., 2012). In theory, the ROI-based RFC measures can be shown to be the sum of the ICA derived RFC both for the within and between networks (Joel et al., 2011;Rosazza et al., 2012).
With ROI-based analysis the brain is first parcellated into predefined anatomical regions, the mean time course for each ROI is then determined. By calculating the temporal correlations in a pairwise fashion between the defined ROIs, for each R-fMRI dataset a correlation coefficient matrix of the ROIs can be obtained for further statistical assessment. Therefore, specific connectivity between specific regions is explicitly tested in a model-driven framework by using the average time courses of the selected ROIs as a temporal model. Since the RFC patterns do not necessarily coincide precisely with the atlas-based ROI definition, all voxels within predefined ROIs are not necessarily a part of the network-of-interest and functionally connected. This can potentially affect the accuracy and sensitivity of the ROI-based analysis (Song et al., 2016). On the other hand, ICA can reveal dynamics and spatially distributed brain networks in a data-driven fashion without the need of a temporal model. Beside motor and sensory networks, ICA studies have identified the brain networks involved in attentional control (Lawrence et al., 2003), including the task-dependent (Vossel et al., 2014) dorsal and ventral lateral attention networks (DAN and VAN), the task-independent (Raichle et al., 2001;Buckner et al., 2008) default mode network (DMN), and the salience network (SN), which was postulated to be the switching control network for the up-regulation of attention networks and the downregulation of the DMN (Menon and Uddin, 2010). The dynamic interactions between the DAN, VAN, DMN, and SN networks are believed to be the key for understanding the function and dysfunction efficient attention allocation for task performance.
Despite the growing consensus regarding the ICA-derived intrinsic RFC networks in the healthy brain with stable spatial components reproduced across studies (Damoiseaux et al., 2006;Smith et al., 2009;Allen et al., 2011), the precise number of independent components (NIC), as a prerequisite input parameter for ICA, is not known a priori. NIC can substantially influence the ICA outcomes (Wang and Li, 2015). Moreover, there is lack of gold standard for the selection of meaningful components to exclude non-interesting noise resources, such as ventricular, vascular, susceptibility, or motion-related artifacts (Wang and Li, 2013).
In this study we refined further of our quantitative data-driven analysis (QDA) framework based on the time course of individual voxel inside the brain. The QDA approach is data-driven as ICA and can generate two types of quantitative RFC metrics for each voxel inside the brain without the need for specifying a particular threshold, model or mode. Since it uses the time course of each voxel within the brain as the reference seed in turn to compute voxel-wise whole-brain correlational coefficient matrix, the size of the correlation matrix is equal to the number of voxels inside the brain. It is typical N > 10 4 for whole-brain R-fMRI datasets with 4 mm voxel size. To facilitate further statistical assessment of the whole-brain correlation matrix, we derive two types of voxel-wise RFC metrics from the correlation matrix, namely the connectivity strength index (CSI) and connectivity density index (CDI). CSI and CDI provide general connectivity metrics of strength and density for the local voxel with the rest of brain, respectively. These metrics can be used for straightforward statistical comparison to assess differences between groups and longitudinal changes of individuals. This is a basic requirement for radiological diagnosis in clinical practice.
Both ICA and ROI-based approaches have previously been applied to study age-related changes in RFC (Bluhm et al., 2008;Biswal et al., 2010;Weissman-Fogel et al., 2010;Zuo et al., 2010;Dennis and Thompson, 2014;Alarcon et al., 2015;Zhang C. et al., 2016). Numerous studies have confirmed that reduced RFC in healthy aging in the DMN is correlated with cognitive deficit (Damoiseaux et al., 2008;Biswal et al., 2010;Campbell et al., 2013;Ferreira and Busatto, 2013;Scheinost et al., 2015). There is accumulating evidence to support the notion that elderly adults typically have reduced RFC across most parts of the DMN, particularly in the dorsal medial prefrontal cortex (mPFC) and the ventral and posterior cingulate cortex (PCC; Campbell et al., 2013;Scheinost et al., 2015). However, in the reported literature there is also considerable variability concerning agerelated RFC differences in the limbic and other DMN subsystems. For example, some studies have found age-related RFC reduction in the hippocampal (Damoiseaux et al., 2008;Campbell et al., 2013;Scheinost et al., 2015) and subcortical regions (Ystad et al., 2010), whereas others reported either no significant decline or elevated RFC in some of the specific hippocampal (Pasquini et al., 2015) and DMN regions (Salami et al., 2014;Damoiseaux et al., 2016). The discrepancies in the reported findings among the different R-fMRI studies may reflect not only variability in the sample characteristics, but also diversity in the data processing methods for deriving the different RFC metrics for connectivity of specific pathways.
The main objective of this study is to develop a QDA framework to analyze R-fMRI data and derive quantitative, model-free, and threshold-free RFC metrics, which are optimally sensitive to physiological and pathological changes in the central nervous systems. We used the proposed metrics to assess if and how adult age in healthy subjects influences these RFC metrics.

Participants
A total of 227 volunteers (aged 18-76 years, male/female = 99/128) completed the study and were recruited into the study through the local media advertisement in the Stockholm region. All participants were right-handed, and native Swedish speakers with normal or corrected-to-normal vision. They all reported being free of a history of neurological, psychiatric, and cardiovascular diseases. None of the participants reported any use of psychotropic drugs. Each subject signed informed consent before completing the magnetic resonance imaging (MRI) examination protocol. They were financially compensated for their participation. The regional ethics committee approved the study protocol 2014/1982-31/1, which was conducted in line with the declaration of Helsinki.

Magnetic Resonance Imaging Data Acquisition Protocol
The MRI data acquisition was conducted on a whole-body 3T clinical MRI scanner (Magnetom Trio, Siemens Medical Solutions, Erlangen, Germany) equipped with a 32-channel phased-array receiving head coil. All data was acquired at Karolinska University Hospital, Huddinge, Stockholm, between noon and 5:00 PM. The MRI data acquisition protocol included the following scanning sessions: (1) 3-plane localizer; (2) Conventional clinical MRI scans including 3D T1-weighted MPRAGE, T2 and FLAIR scans; and (3) A session of 375 s long R-fMRI measurements. The main acquisition parameters for the R-fMRI data included the following: TE/TR 35/2,500 ms, flip angle = 90 • , 34 slices of 3.5 mm thick, FOV = 225 mm, matrix size = 76 × 76, data acquisition acceleration with GRAPPA parallel imaging method (iPAT = 2), and 150 dynamic timeframes. The T1-weighted MPRAGE images used for co-registration with functional images were acquired with the following parameters: TR = 1,900 ms, TE = 2.52 ms, FA = 9 degrees, FOV = 256, voxel size 1 mm × 1 mm × 1 mm. The acquisition parameters for the FLAIR image were the following: TE/TR = 89/9,000 ms, flip angle = 130 • ; inversion time (TI) = 2,500 ms, slice thickness = 4.0 mm, and FOV = 199 × 220 mm. An experienced radiologist inspected both the FLAIR and T1-weighted images for potential signs of neuropathology.
We used foam patting to fix each subject's head carefully in the head coil to reduce involuntary head motions. During the R-fMRI data acquisition the participants were instructed to focus their sight on a white cross in black background projected on a screen installed in front of their eyes. The subjects were also instructed to not think about anything particular during the R-fMRI session.

Resting-State Functional Magnetic Resonance Imaging Data Pre-processing
The R-fMRI datasets underwent a preprocessing procedure, which has been described elsewhere in details (Li et al., 2021) and was performed with AFNI (Version Debian-16.2.07∼dfsg.1-3∼nd14.04+1, http://afni.nimh.nih.gov/afni) and FSL 1 programs with a bash wrapper shell Li, 2013, 2015). After temporal de-spiking, six-parameter rigid body image registration was performed for motion correction. The average volume for each motion-corrected time series was used to generate a brain mask to minimize the inclusion of the extra-cerebral tissues. Spatial normalization to the standard MNI template was performed using a 12-parameter affine transformation and a mutual-information cost function. During the affine transformation the imaging data were also re-sampled to isotropic resolution using a Gaussian kernel with 4 mm full width at half maximum (FWHM). The coregistered average image volume for the cohort has 28,146 non-zero voxels inside the brain and was used to generate the average brain mask for the preprocessed whole-brain R-fMRI data with 4 mm spatial resolution. Nuisance signal removal was performed by voxel-wise regression using 14 regressors based on the motion correction parameters, average signal of the ventricles and their 1st order derivatives. After baseline trend removal up to the third order polynomial, effective band-pass filtering was performed using lowpass filtering at 0.08 Hz. Local Gaussian smoothing up to FWHM = 4 mm was performed using an eroded gray matter mask (Wang and Li, 2015).
Pearson's correlation coefficients (CC) were computed between the time courses of all pairs of voxels inside the brain, leading to a whole-brain functional connectivity matrix for each subject. This computation was performed for all voxels located within the brain mask, which was generated by overlapping the registered brains of all participants. This brain mask contained 28,146 voxels and each voxel inside the brain was used as the seed voxel in turn. Therefore, the size of the CC matrix size is 28,146 × 28,146. Each row or column of the CC matrix corresponds to the CC image volume for the seed voxel with the rest of the brain. That is the connectivity map for the seed voxel. As schematically illustrated in Figure 1, based on the CC histogram for each row of the matrix we derived the following two types of FIGURE 1 | A schematic overview to illustrate the QDA framework. With QDA the time course of each voxel is used in turn to compute the whole-brain CC matrix. For each row of the CC matrix, we compute a CC histogram with 200 evenly binned intervals within [−1, 1]. The histogram shown in the graph is the cohort's average CC histogram for a voxel within the PCC as marked with the cross. Two types of RFC images are derived from the CC matrix: (1) CSI P and CSI N whose voxel values are the averages of the positives and negatives in each row of the CC matrix, respectively. (2) CDI P and CDI N whose voxel values are the positive and negative parts of the convolution between the CC histogram and the kernel, respectively. threshold-free voxel-wise RFC metrics: the CSI and CDI. As we are interested in systematically investigating all relevant synchronized activities in the whole brain, we quantify the negative and positive portions of the CC histogram separately to avoid information cancelation, sensitivity reduction, and statistical interference. From here on, the subscripts "N" and "P" are used to indicate the negative and positive portions of the RFC metrics, respectively. The metrics without subscripts refer to the mixed measures without distinction of the negative and position correlations.
As shown in Figure 1, the voxel value for the CSI P and CSI N are defined as the averages of the positives and negatives in each row of the CC matrix, respectively. That is Where CC row refers to a row in the CC matrix. np and nn refer to the number of positive and negative correlation coefficients in a row of the CC matrix, respectively. The voxel values for CDI are defined as the convolution between the CC histogram and a kernel function. That is The CDI P and CDI N correspond to the positive and negative portions of the convolution defined in Eq. (3), respectively. To facilitate statistical comparison it is useful to transform the raw RFC metrics into standard Z-score using the following formula: Where µ and σ are the mean and standard deviation of the corresponding RFC metrics, respectively. For optimization of the CDI sensitivity, we investigated 6 different kernel functions including where x⊂ [−1,1] corresponds to the interval of the correlation coefficients. The kernels are also graphically depicted in Figure 2.
Frontiers in Neuroscience | www.frontiersin.org FIGURE 2 | The six different kernel functions investigated in the study to derive the CDI P and CDI N metrics. The widely used threshold method can be considered as the case for the square-well kernel function (k 6 ).
The kernel should weight the higher correlation coefficients more than the lower ones. The widely used threshold approach can be considered as the case of the square-well kernel function k 6 . For illustration, an arbitrary threshold of 0.3 was used here. The CSI metrics can also be considered as a special case of CDI corresponding to a kernel of the sign function.

Statistical Analyses
To investigate if and how the RFC metrics are influenced by heathy aging for the studied cohort we performed voxel-wise linear regression analyses of the CSI and CDI metrics versus the subject's age, while gender was treated as a covariate by using the AFNI program 3dRegAna to extract the regression parameter β and linear coefficient r. The statistical significance was assessed by using a two-step approach. Firstly, we imposed a voxel-wise threshold p < 0.001 (uncorrected corresponding t-score ≥ 3.34) to form the initial cluster candidates. Secondly, we performed permutation simulations without assuming a particular form of probability distribution for the voxel values in the statistic images to identify the brain ROI out of the initially detected clusters at family-wise error rate p ≤ 0.05. Using the detected ROIs as masks, the mean values of the RFC metrics for each ROI were evaluated and plotted against the subjects' age. Besides linear regression analysis with age, we performed also verification using two-sample t-test between the young and elderly subgroups. For this, we selected all subjects aged 18-30 years as the young subgroup (n = 124, males/females = 51/73), and all subjects aged 64-76 years as the elderly subgroup (n = 76, males/females = 35/41). To keep sufficient age gap between the young and elderly subgroups the remaining 27 subjects in the age range of 31-63 years old were excluded from the t-test. In the selection of subgroups we attempted to minimize the number of excluded subjects with intermediate ages, maximize the age gap between the subgroups, and keep similar number of subjects and age ranges. It should be emphasized that all 227 subjects were included in the regression analysis.

The Quantitative Data-Driven Analysis Framework
The CC histogram for each seed voxel in the brain is dependent on its location in the brain (see Supplementary Materials). Figure 3 shows the average CC histogram of the cohort for a seed voxel in the PCC as illustrated by the cross in Figure 1.
The histogram is somewhat asymmetric and shifted toward the positive side. This is quite typical at least for voxels within gray matter. Selecting different threshold values along the histogram allows us to examine the RFC networks of different connection strengths associated with the selected seed voxel. As shown in Figure 3, at high negative threshold (Figures 3A,B) we observe the DMN. At low negative threshold, we observe its association with cerebral spinal fluid (CSF) space and white matter (Figure 3D). At moderately high positive threshold, the PCC voxel is not only a part of the DMN, but also connected to most of the cortical gray matter (Figures 3E,F). At high positive threshold, the PCC voxel is associated with the posterior portion of the DMN and the visual cortex ( Figure 3C). The visual cortex looks relatively bright at high threshold indicating a relatively high number of voxels are associated with the visual network or voxels within the visual cortex are associated with each other at high threshold criterion. The idea of the QDA framework is to avoid the arbitrary threshold and optimize the contribution of meaningful informatics to the quantitative RFC metrics.

Figure 4
shows an axial slice of CDI P and CDI N images for a typical R-fMRI dataset (from a 36 year old male subject). Multiple brain regions depict high CDI P including the bilateral mPFC, superior and middle temporal gyri (MTG), inferior and superior parietal lobule, precuneus and PCC. These regions have been described as RFC hubs implying their important role in neural signaling and communication across the brain (Buckner et al., 2009;Tomasi and Volkow, 2011). On the other hand, the PCC, insula cortex. White matter and CSF regions have high CDI N metric. The contrast and intensity variations across each row in Figure 4 demonstrate that selection of the kernel FIGURE 4 | An axial slice of the CDI P (upper row) and CDI N (lower row) metrics derived from a typical R-fMRI dataset (a male subject of 36 year's age). The images from left to right depict the results for the following 6 kernel functions | x|, | x 2 |, | x 3 |, | x 4 |, sin 2 (π/2x), and step(| x| −0.3), respectively. function can optimize the contrast and signal-to-noise ratio of the CDI metrics.

Resting-State Functional Connectivity Changes Associated With the Adult Age
The linear regression results for the CSI, CSI N , and CSI P data versus subjects' age are summarized in Figure 5 and Table 1. The corresponding results for the CDI, CDI N , and CDI P are shown in Figure 6 and Table 2. The CSI metric without separation FIGURE 5 | Brain regions with significant correlation (p < 0.05, corrected) between the connectivity strength metrics and the subject's age. The results for the CSI (A), CSI N (B), and CSI P (C) are depicted separately. The Color bar shows the t-score level.
of the negative and positive correlations shows decline of the functional connectivity strength with age in the superior and middle prefrontal gyrus (MFG) and increase of connectivity strength in the precuneus and right inferior parietal lobule (r-IPL). The more specifically defined CSI N and CSI P metrics are more sensitive to the adult age effect and the detected brain volumes with significant aging effect are nearly tripled compared with that for the CSI metric. With CSI N and CSI P we also observe a more intricate pattern of change with the adult age, which are summarized as follows: (1) The CSI P shows mainly decline trend with adult age (negative β and r) in the extended DMN including superior and MFG, PCC, bilateral insula cortex and left middle temporal gyrus (l-MTG) except for putamen where upregulation of CSI P was observed.
(2) The CSI N depicts a more complicated pattern of dependence on adult age. The negative connectivity strength was reduced (positive β and r) with adult age in the PCC, right insula cortex and IPL, while enhancement (negative β and r) was detected in the sensorimotor network (paracentral lobule, bilateral postcentral gyri), bilateral parahippocampal cortices (PHC), and right superior temporal gyrus. (3) There are two brain regions where both the CSI N and CSI P demonstrated significant reduction trend with adult age, which were detected by applying the logical "AND" operation to the regression results for the CSI P and CSI N . As shown in Table 2 and Figure 7, the two overlapping ROIs in the PCC and r-insula cortex depict significant down-regulation of CSI P and CSI N metrics with the subjects' age.
To study the specific connectivity associated with the two ROIs defined by the overlap between the CSI N and CSI P metrics, we computed Pearson's correlation maps for the time courses of the seeds as defined by the overlapping ROIs depicted in Figure 7A. As expected, the associated RFC network for the PCC ROI is obviously the well-known DMN and include The volume, center of mass coordinates in MNI space, regression parameter (β), linear correlation coefficient (r), statistical significance (p), and anatomic annotations are specified. The default is bilateral, while R and L-indicate the right and left hemisphere of the brain, respectively. CSINP indicates the overlapping results between CSIN and CSIP.
4 negatively correlated brain regions, which are the bilateral IPL and insula cortices (see Figure 7B). On the other hand, the associated RFC network for the insula ROI includes the PCC and bilateral precuneus as the negatively correlated brain regions ( Figure 7C) Figure 8. Shows the anti-correlated brain regions between the above 2 RFC networks as obtained by multiplying the two correlation maps with each other and applying negative threshold at CC ≤ (−0.5). The mutually inclusive anti-correlation between the PCC and the right insular cortex explains why both CSI P and CSI N metrics in these regions depict declines with the adult age. Figure 9 shows the ROI average of the CSI N and CSI P metrics in the PCC and right insula cortex as a function of the subject's age. With normal aging, both the CSI P and CSI N are reduced in these brain regions (overlap shown in Figure 7A). Therefore, the PCC and right insula are particularly sensitive to the adult age effect. However, the aging effect is barely detectable by the unseparated CSI metric (see Table 1).
As expected, the CDI P and CDI N metrics derived by using the different kernels differ in their sensitivity in detecting the adult age effect. Figure 10A shows the detected brain volumes where the CDI P and CDI N metrics are significantly associated with the adult age. The sensitivity difference of the kernels is also manifested in the regression parameter β which are detailed in Table 3 and Figure 10B. To compare the similarity of the detected aging effects among the CDI metrics of different kernels, we assessed the joint overlapping brain regions detected by the different CDI N and CDI P metrics of different kernels. The observed overall trends of RFC enhancement or decline with age are quite similar. The joint overlapping volumes for the CDI P and CDI N metrics of different kernels are 733 and 671 voxels, respectively. Moreover, there is also a reasonable anatomic consistency between the results of the connectivity strength metrics and connectivity density metrics. As detailed in Tables 2, 3, the anatomical locations of the joint overlapping regions for the different CDI P metrics match those for the 3 largest ROIs identified by the CSI P results (see Table 1). Similarly, the brain regions of the joint overlapping for the different CDI N metrics are largely the same as those identified by the CSI N data (see Table 1). However, it should be noted that the β parameters for the CDI N and CSI N have opposite signs even through the trend of change with the adult age is the same. This is because the negative connectivity strength (CSI N ) is negative in nature, while the connectivity density corresponding to the negative correlation (CDI N ) is always positive. Therefore, the enhancement of the negative connectivity strength (CSI N ) with age (for example in the sensorimotor network) corresponds to a negative β, while the connectivity density result corresponds to a positive β value.

Effects of Adult Age on Resting-State Functional Connectivity
Age is an important risk factor for declines of neural cognitive functions and pathology of neurodegenerative diseases. It is also FIGURE 6 | Brain regions with significant correlation (p < 0.05, corrected) between the connectivity strength metrics and the subject's age. The results for the CDI (A), CDI N (B), and CDI P (C) are depicted separately. The Color bar shows the t-score level. a complex metric difficult to precisely interpret the involved physiology. Healthy individuals of similar age may have quite different vascular and brain-health status. It follows that age is not a single strongest predictor for the RFC in the brain. This is likely to be the reason why the linear regressions of the RFC metrics with the adult age depict substantial scatters and relative low correlation coefficients. The impact of the potential confounds and pre-processing strategies that can mitigate them have been extensively investigated in the published literature (Bluhm et al., 2008;Biswal et al., 2010;Weissman-Fogel et al., 2010;Zuo et al., 2010;Dennis and Thompson, 2014;Alarcon et al., 2015;Zhang C. et al., 2016;Geerligs et al., 2017;Hussein et al., 2020). Here we focus on comparing our findings in the context of documented literature results, particularly the adult age effect in the DMN, dorsal attention network (DAN), sensorimotor network, and subcortical brain regions.
With QDA, we found support for RFC decline with advancing adult age in multiple brain regions of the DMN and DAN, including superior and MFG, PCC, MTG, and IPL. Age-related RFC decrements in the DMN and DAN have previously been reported in numerous R-fMRI studies using ROI and ICA based analysis (Buckner et al., 2008;Damoiseaux et al., 2008;Dennis and Thompson, 2014;Scheinost et al., 2015;Luo et al., 2020). Our findings regarding to the RFC changes in the DMN are overall in agreement with previous reported results (Damoiseaux et al., 2008(Damoiseaux et al., , 2016Koch et al., 2010;Ystad et al., 2010;Williams, 2013;Lu et al., 2014;Persson et al., 2014Persson et al., , 2015Salami et al., 2014). Besides the DMN and DAN, normal aging was associated with RFC increase in the sensorimotor, subcortical network, and para-hippocampal cortex. This has also been reported previously (Persson et al., 2015(Persson et al., , 2016Damoiseaux et al., 2016;Geerligs et al., 2017;Hussein et al., 2020;Luo et al., 2020). We didn't find significant age-related RFC declines in precuneus and specific sub-regions of the hippocampal cortex as reported in previous studies (Salami et al., 2014;Damoiseaux et al., 2016). Since we assessed the negative and positive correlation separately, this may allow us to detect more intricate agerelated RFC changes in the brain. To illustrate this point, we analyzed further the 3 ROIs with significant correlation between the CSI and the subject's age. As shown in Tables 1, 4 and Figure 11, the detected ROI in the precuneus depicted significant positive linear correlation between CSI and the subject's age (β = 9.50 × 10 −3 , r = 0.459), even though the CSI P and CSI N in the same ROI showed only a slight (not significant) increment and decrement with age, respectively, i.e., contribution from a low-significant CSI P increment and a non-significant CSI N decrement resulted in a highly significant increment trend in the CSI metric. With the same line of reasoning, we can explain why the MFG ROI detected by the CSI metric is much smaller than that detected by the CSI P metric, because the decremental trend in the CSI P metric was partially canceled by the CSI N contribution. This can also explain why we didn't detect significant CSI decrement with the adult age in the PCC and R-insula, because both the CSI P and CSI N metrics exhibited significant decremental trends with age and their contributions annulled each other. Therefore, it is important to pay attention to the precise definition of the RFC when comparing the results of different studies.
Both CSI P (r = 0.506, see Table 1) and CDI P (r = 0.577, see Table 2) showed age-related enhancement Caudate/putamen and the association are quite strong. This finding based on QDA approach are consistent with previous reports (Manza et al., 2015;Rieckmann et al., 2018) from ROI-based studies aimed to investigate the aging effect on specific functional connectivity of in the striatum-cortical system. It is well known that the striatum-somatomotor connection is primarily associated with motor performance, especially the "automatic" performance of already learned movements. It has been reported that posterior putamen and pallidum decreases in connectivity to left somatomotor cortex with age (Manza et al., 2015). This provides a reasonable explanation for the commonly observed motor deficits as in elderly subjects. There is also growing evidence support the notion that striatum-cortical connectivity is potentially important for memory function at older age. Intriguingly, this is related to the enhanced RFC with age in the putamen/caudate. Several studies have reported that increased striatum functional connectivity in older adults typically reflects less negative connectivity between two regions belonging to different networks, and the increased connectivity is often negatively associated with cognitive performance (Rieckmann et al., 2018). Therefore, better understanding the RFC change with The volume, center of mass coordinates in MNI space, regression parameter (β), linear correlation coefficient (r), statistical significance (p), and anatomic annotations are specified. The default is bilateral, while R and L-indicate the right and left hemisphere of the brain, respectively. CDI NP indicates the overlapping results between CDI N and CDI P .
adult age in the striatum-cortical system can be potentially useful for assessing motor as well cognitive functions in elderly subjects.

Methodological Issues
The QDA framework proposed in the study is a voxel-wise and data-driven approach. It has the following two unique features: (1) It can avoid confounding caused by the cancelation of the negative and positive correlations by assessing the negative and positive portions of the CC histogram separately; (2) It derives different RFC metrics based on the connectivity strength and density by utilizing the concept of convolutions with different kernels. The metrics weight all the correlations of a given voxel with the rest of the brain according to the amplitudes of the correlation coefficients and disregard the anatomical distance between the correlation pairs. This permits a comprehensive characterization of the intrinsic activities of each voxel without the use of an arbitrary threshold. The QDA approach can encapsulate the widely used threshold approach as a special case of the square-well kernel function. The widely used degree centrality corresponds precisely this square-well kernel situation which adopts a somewhat arbitrary threshold and every connection above the threshold are weighted equally. Even the CSI metrics can be encapsulated under the convolution concept for a special kernel of the sign function. This provides a unified view for RFC and can facilitate its further optimization. The QDA framework uses the time course of each voxel within the brain as the seed reference to compute voxel-wise whole-brain correlational coefficient matrix. For whole-brain R-fMRI data acquired at 4 mm spatial resolution, the correlational coefficient matrix is in the order of 10 5 and is currently not practical for direct visualization and statistic assessment. Particularly, when data are acquired with higher spatial resolution, e.g., 2 mm, the matrix size is increased by 8 × 8 times. Therefore, for data reduction, we derived two types of voxel-wise RFC metrics from the correlational coefficient matrix without the need for specifying a particular seed, threshold, model, or mode. As their names indicated, CSI and CDI are aimed to capture the local (voxel) connectivity strength and density with rest of brain, respectively. The QDA metrics can assess the general connectivity with the rest of the brain without specifying a specific path or network. The QDA method does not highlight the specific FIGURE 7 | The overlapping ROIs in the PCC and right insula cortex where both the CSI P and CSI N metrics depict significant decline with the adult age (A). The one sample t-test maps for the Pearson's correlation maps associated the seeds defined as the overlapping ROI in the PCC (B) and insula cortex (C).
connectivity changes between selected brain regions. The precise neural correlate of R-fMRI signal is currently not well understood (Hyder and Rothman, 2010). However, it is reasonable to assume that the R-fMRI signal fluctuations indirectly reflect the slow modulations of neuronal activities at rest. Furthermore, the sigmoid function has been widely accepted as the logistic function of neuronal activation instead of a square-well. With current status of knowledge, we cannot identify a convolution kernel to reveal a particular feature of the neurophysiology. However, we can attempt to optimize the kernel to reduce bias and improve sensitivity of the RFC metrics to pathophysiological changes. The current results based on the QDA framework should be interpreted in the context of some technical and biological limitations. Firstly, at a TR of 2,500 ms, the cardiac and respiratory fluctuation effects might be aliased into the low frequency R-fMR signal fluctuations. The regression up-to the 1st order derivative of the head motions and lowpass filtering could not eliminate the effects of these physiological noises (Muschelli et al., 2014;Pruim et al., 2015;Bright et al., 2017;Parkes et al., 2018). Thus, these aliasing effects could reduce the specificity of the RFC metrics, or even might further confound the detected RFC differences between the young and elderly subgroups. With the more up-to-date acquisition techniques, such as multi-band simultaneous acquisition of multiple slices and compress-sensing with high under-sampling factor, it is possible to use a shorter TR (e.g., 500 ms) and higher spatial resolution for the data acquisition. Therefore, these physiological effects may be further mitigated.
Secondly, the resting state is associated with spontaneous thoughts and cognitive processing, we cannot exclude the possibility that differences in spontaneous thoughts may exist between the young and elderly subjects (Wu et al., 2007). However, considering the overall consistency of our results with the previous studies, particularly the results from the longitudinal studies (Fjell et al., 2014;Ng et al., 2016;Staffaroni et al., 2018;Li Q. et al., 2020), it is unlikely that these differences have major FIGURE 9 | Scattered plots of the regression against age for the overlapping results between the CSI N and the CSI P metrics (see the bottom rows of Table 1). The ROI average of the CSI N metric against the subject's age for the overlapping ROI in the PCC (A). The ROI average for the CSI N metric against the subject's age for the overlapping ROI in the right insula cortex (B). The ROI average for the CSI P metric against the subject's age for the overlapping ROI in the PCC (C). The ROI average for the CSI P metric against the subject's age for the overlapping ROI in the insula cortex (D). The lines show the linear regression results of the RFC metrics against the subject's ages.
influence on our findings. These initial findings encourage the future use of QDA as a tool to analyze longitudinal R-fMRI data aimed to develop a comprehensive understanding of age-or pathology-related brain functional changes.
Thirdly, the generalizability, or external validity issue should be considered. This is due to the non-random recruitment procedures and relying on a sample of convenience. The sample size used in this study (N = 227) is moderate, includes unbalanced young and elderly subgroups reflecting the difficulties to recruit elderly healthy subjects. The ages of the participants range from young to old adulthood (reflecting the age of participants in most neuroimaging studies). The age-related RFC differences observed in this study were relatively small but quite robust. However, the results from this cross-sectional study of the cohort cannot distinguish whether the RFC changes in the brain regions are due to gradual changes throughout the adulthood or a more sudden change at later stage in life.
Fourthly, the R-fMRI data were acquired under open-eye condition. Recent studies indicate that opening versus closing eyes at resting-state results in RFC difference between V1 with DMN and SNs (Costumero et al., 2020). This may explain why we did not detect significant RFC change with age in the visual cortex.

Negative Cross Correlation in White Matter and Cerebral Spinal Fluid
As discussed above negative correlation is an important fraction of the CC histogram irrespective of the tissue type and anatomical location of the voxel in question. In the published literature, there is also a rapid growing interest in studying the negative correlations between the voxels Weissenbacher et al., 2009;Bianciardi et al., 2011;Schwarz and McGonigle, 2011;Gonzalez-Castillo et al., 2012;Gopinath et al., 2015;Liu et al., 2015;Spreng et al., 2016;Chen et al., 2020). It is clear that the negative portion of the CC histogram is more dominant for voxels in CSF (Gruszecki et al., 2018) and white matter (McColgan et al., 2017;Gore et al., 2019; FIGURE 10 | The total volumes of the detected brain regions with significant correlation (p < 0.05, corrected) between the connectivity density index (CDI P and CDI N ) and the subject's age as a function of the kernels (A). The average regression parameter β for the detected brain regions as a function of the kernels (B). The negative and positive correlations were assessed separately. The volume, center of mass coordinates in MNI space, and anatomic annotations, regression parameters (β), linear correlation coefficient (r), statistical significance (p), and anatomic annotations are specified. The default is bilateral, while R and L-indicate the right and left hemisphere of the brain, respectively. CDI N and CDI P indicate the joint overlaps among the CDI N and CDI P metrics of the different kernels, respectively. The β, r, and p are the average results for the 6 different kernels. Wu et al., 2019;Li M. et al., 2020). However, the negative portion cannot be ignored even for voxels in the gray matter (see Supplementary Materials). To avoid confound caused by inappropriate preprocessing pipelines, we have carefully tested and updated our preprocessing pipeline. We did not implement the global signal regression (GSR) which removes the mean signal averaged over the entire brain. GSR removal via linear regression is one of the most controversial procedures in the analysis of R-fMRI data Weissenbacher et al., 2009). On one hand, the global mean signal contains variance associated with respiratory, scanner-, and motion-related artifacts. Its removal by GSR can improve various quality control metrics, which enhances the anatomical specificity of RFC networks, and increase the explained behavioral variance. On the other hand, GSR alters the distribution of regional signal correlations in the brain, can induce artefactual anti-correlation patterns, may remove real The CSI P and CSI N results are based on the masks determined solely by the CSI results.
FIGURE 11 | The ROI average of the CSI P , CSI N , and CSI metrics against the subject's age for the 3 ROIs with significant correlation between CSI and the subject's age. The details of the regression results are summarized in Table 4. The columns 1 to 3 are the results for the ROIs in the precuneus, MFG, and R-IPL, respectively. The rows 1 to 3 are the results for the CSI P , CSI N and CSI metrics, respectively. The ROI masks are solely based on the CSI metric only.
neural signal, and can distort RFC metrics. The brain masked "global signal" is usually misunderstood, because it is not "global" and its variance contains dominant contributions from different domains of the voxels with temporally coherent signal variation.
To limit the study in a reasonable scope, in the discussion of the adult age effect on RFC we focused on gray matter and did not discuss white matter and CSF related issues. However, it should be pointed out that aging effects in white matter (McColgan et al., 2017;Gore et al., 2019;Wu et al., 2019;Li M. et al., 2020) and CSF (Sakka et al., 2011;Gruszecki et al., 2018) are also worth exploring. There is indeed a rapid growing interest in these arenas in published literature (Sakka et al., 2011;McColgan et al., 2017;Gruszecki et al., 2018;Gore et al., 2019;Wu et al., 2019;Li M. et al., 2020), particularly in the context of the age effect for the glymphatic system.

CONCLUSION
The proposed QDA framework can provide data-driven, threshold-free and voxel-wise analysis of R-fMRI data and offer a unified view for RFC metrics which can facilitate further development and optimization of the RFC metrics by choosing appropriate kernel functions. The QDA results for the adult age effect are largely consistent with previously published results based on other analysis methods. Moreover, our new findings based on the separate assessment of the negative and positive correlations can improve the sensitivity of the RFC metrics to physiological changes associated with the advancing adult age and may clarify some of the confounding reports in the literature regarding to the DMN and sensorimotor network involvement in normal aging.

DATA AVAILABILITY STATEMENT
The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: http://dx.doi.org/10.17 632/pt9d2rdv46.1.

ETHICS STATEMENT
The studies involving human participants were reviewed and approved by the regional Ethics Committee in Stockholm and the study protocols included 2012/1511-31/2 and 2014/1982-31/1. The patients/participants provided their written informed consent to participate in this study.